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Abstract 

We measure three-particle Bose-Einstein correlations in hadronic Z decay with 
the L3 detector at LEP. Genuine three-particle Bose-Einstein correlations are ob- 
served. By comparing two- and three-particle correlations we find that the data are 
consistent with fully incoherent pion production. 
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Introduction 



So far, no theory exists which can describe the non-perturbative process of hadron production 
in general and Bose-Einstein (BE) effects in particular. The latter are expected from general 
spin statistics considerations. To help understand these phenomena, studies of identical-boson 
correlations in e + e~ collisions at LEP have been performed in terms of the absolute four- 
momentum difference Q []]. as well as in two- and three-dimensional distributions in components 

of q mm. 

It has long been realized that the shape and size in space-time of a source of pions can be 
determined, as a consequence of the interference of identical bosons, from the shape and size 
of the correlation function of two identical pions in energy-momentum space j3j. Additional 
information can be derived from higher-order correlations. Furthermore, such correlations 
constitute an important theoretical issue for the understanding of Bose-Einstein correlations 
(BEC) 0- 

In this Letter three-particle correlations are analysed. These correlations are sensitive to 
asymmetries in the particle production mechanism PIT] which cannot be studied by two-particle 
correlations. In addition, the combination of two- and three-particle correlation analyses gives 
access to the degree of coherence of pion production [HUH!, which is very difficult to investigate 
from two-particle correlations alone due to the effect of long-lived resonances on the correlation 
function. The DELPHI ^Uj and OPAL [TT] collaborations have both studied three-particle 
correlations but did not investigate the degree of coherence. 

The Data and Monte Carlo 

The data used in this analysis were collected by the L3 detector ^2] i n 1994 at a centre-of- 
mass energy of 91.2 GeV and correspond to a total integrated luminosity of 48.1 pb _1 . The 
Monte Carlo (MC) event generators JETSET ^Hj and HERWIG ^1] are used to simulate the 
signal process. Within JETSET, BEC are simulated using the BE algorithm H3Q6 The 
generated events are passed through the L3 detector simulation program, which is based on 
the GEANT [T7] and GHEISHA [TH] programs, reconstructed and subjected to the same 
selection criteria as the data. 

The event selection is identical to that presented in Reference El resulting in about one 
million hadronic Z decay events, with an average track multiplicity of about 12. Two additional 
cuts are performed in order to reduce the dependence of the detector correction on the MC 
model used: tracks with measured momentum greater than 1 GeV are rejected, as are pairs of 
like-sign tracks with opening angle below 3°. This results in an average track multiplicity of 
about 7. For the computation of three-particle correlations, each possible triplet of like-sign 
tracks is used to compute the variable Q 3 = \JQ\2 + Q23 + Qiii where = \J —(pi — Pj) 2 is 
the absolute four- momentum difference between particles i and j. Since Qtj, and thus Q 3 , 
depends both on the energy of the particles and on the angle between them, small Qij can 
be due to small angles or low energies. In a MC generator with BE effects, the fraction of 
pairs at small with small angle is larger than in one without. Consequently, the estimated 
detection efficiency depends on the MC model used. The momentum and opening angle cuts 
reduce this model dependence. After selection, the average triplet multiplicity is about 6. In 

The Bose-Einstein simulation is done by the subroutine LUBOEI, with the values PARJ(92) = 1.5 and 
PARJ(93)=0.33GeV. 
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the region of interest, Q 3 < 1 GeV, the loss of triplets by the momentum and opening angle 
cut is about 40%. 

The momentum cut improves the resolution of Q 3 by a factor three with respect to that 
for the full momentum spectrum. Using MC events, its average is estimated to be 26 MeV for 
triplets of tracks with Q 3 < 0.8 GeV. We choose a bin size of 40 MeV, somewhat larger than 
this resolution. 

In Figure ^ the data are compared to JETSET (with and without BE effects) and HERWIG 
(not having a BE option) at the detector level, after performing all the cuts mentioned above, 
in the three-particle distributions J28<p, Y^&Q-i an d Qz- The sums run over the three pairs 
of like-sign tracks in the triplet and 8(f) and 89 are the absolute differences in azimuthal and 
polar angle between two tracks, respectively. Within 10%, the angular distributions of the 
MC models agree with those of the data. None of the models describes the Q 3 distribution: 
JETSET with BE effects overestimates the data by approximately 20% at low Q 3 , even though 
we found good agreement for 8(f), 86, and Q |2J. JETSET without BE effects and HERWIG 
grossly underestimate the data at low Q 3 . The statistics for Q 3 < 160 MeV are so poor, that 
this region is rejected from the analysis. 



The Analysis 

The three-particle number density P3(pi,p 2 ,P3) of particles with four-momenta Px,P2 an d P3 
can be described in terms of single-particle, two-particle and genuine three-particle densities as 

P3(P1,P2,P3) = Pl(Pl)pl(P2)pl(p3)+J2 (Pl(Pl) [P2(P2,P3) ~ Pi (^)pi (ps)} } + C 3 (pi , P2, P3) , (1) 

(3) 

where the sum is over the three possible permutations and C 3 is the third-order cumulant, which 
measures the genuine three-particle correlations. The p\p 2 terms contain all the two-particle 
correlations. In order to focus on the correlation due to BE interference, we replace products 
of single-particle densities by the corresponding two- or three-particle density, po, which would 
occur in the absence of BEC, and define the correlation functions 

„ f v _ f>2(pi,P2) „ ( P3{PI,P2,P3) ,„x 

R2{PUP2 ) = —( r , R3{PUP2,P3 ) = —( 7 • 2 

P0{P1,P2> Po{Pl,P2,P3> 

Assuming the absence of two-particle correlations, i.e., P2(pi,P2) = Pi (pi)Pi(j>2), results in 

Rf^ e ( Pl ,P2,P3) = 1 + C3(P1 ' P2 ' P3) . (3) 

P0{Pl,P2,P3) 

The kinematical variable Q 3 is used to study three-particle correlations. For a three-pion 



system, Q 3 = \J Mf 23 — 9m^, with M 12 3 the invariant mass of the pion triplet and m n the mass 
of the pion. In this Letter, p 3 is defined as 

with N cv the number of selected events and n tr ipiets the number of triplets of like-sign tracks, 
and p2 is defined analoguously. 
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Assuming totally incoherent production of particles and a source density f(x) in space-time 
with no dependence on the four-momentum of the emitted particle, the BE correlation functions 
is related to the source density by [HUH] 



Rf nmne (Q 12 , Q 23 , Q31) 




1 + 

1 + |F(Q 12 )| 2 + |F(Q 23 )| 2 +|F(Q 31 )| 2 

+ 2Re{F(Q 12 )F(Q 23 )F(Q 31 )} 
l + 2Re{F(Q 12 )F(Q 23 )F(Q 31 )} 



(6) 
(7) 



(5) 



where F(Qij) is the Fourier transform of f(x). 

R2 does not depend on the phase <pij contained in F{Qij) = \F(Qij)\ exp(i0jj). However, 
this phase survives in the three-particle BE correlation functions, Equations © and ((7j). As- 
suming fully incoherent particle production, the phase <pij can be non-zero only if the space-time 
distribution of the source is asymmetric and Qij > 0. Defining 



then for an incoherent source Equations © and (|IJ) imply that u = cos0, where = 0i 2 + 
023 + 03i- Furthermore, as Qij — > 0, then 0^ — > 0, and hence u> — > 1. For Qij > 0, a deviation 
from unity can be caused by an asymmetry in the production. However, this will only result 
in a small (a few percent) reduction of u jHllZl, and this only in the case where the asymmetry 
occurs around the point of highest emissivity. It is important to emphasize that for (partially) 
coherent sources, u can still be defined by Equation (JHJ), but Equations (jSHZI) are no longer 
valid, in which case more complicated expressions are needed [Zj , and one can no longer deduce 
that u = cos0 or that u — > 1 as Qij — > 0. In at least one type of model, one can make the 
stronger statement that the limit u — 1 at Qij — > can only be reached if the source is fully 
incoherent [20 . 



The reference sample, from which p is determined, is formed by mixing particles from different 
data events in the following way Firstly, 1000 events are rotated to a system with the z-axis 
along the thrust axis and are stored in a "pool" . Then, tracks of each new event outside the pool 
are exchanged with tracks of the same charge from events in the pool having about the same 
(within about 20%) multiplicity, under the condition that all tracks originate from different 
events. Thus, after this procedure the new event consists of tracks originating from different 
events in the pool, and its original tracks have entered the pool. This updating process prevents 
any regularities in the reference sample. Finally, Q 3 is calculated for each triplet of like-sign 
tracks, resulting in the density p m i x . 

This mixing procedure removes more correlations than just those of BE, e.g., those from 
energy-momentum conservation and from resonances. This effect is estimated using a MC 
model with no BE effects (JETSET or HERWIG) at generator level and using pions only. 
Thus, in the absence of BEC, the corrected three-particle density is given by 



^(Ql2, Q23, Q31) 



Rf nmne (Q 'n, Q23, Q31) — 1 



(8) 



2yJ(Ba(Qi2) ~ WMQw) 1 l)(^ 2 (Q 3 i) - 1) 



Determination of R3 and R3 



Po(<?3) = Pmi x (<23)C mi x(<23), where C mix (<2 3 ) 



p 3 (Q 3 ) 

Pmix(£t?3) 



J MC,noBE 



(9) 
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The density p 3 , measured in data, must be corrected for detector resolution, acceptance, 
efficiency and for particle misidentification. For this we use a multiplicative factor, det, derived 
from MC studies. Since no hadrons are identified in the analysis, Cdet is given by the ratio of the 
three-pion correlation function found from MC events at generator level to the three-particle 
correlation function found using all particles after full detector simulation, reconstruction and 
selection. Combining this correction factor with Equations (j2j) and (jHJ) results in 

The genuine three-particle BE correlation function, R& enmne ^ j s obtained via 

Rf~ = R 3 -R h2 + l , (11) 

where Ri 2 = P1P2) / Po — 2 is the contribution due to two-particle correlations, as may be 
seen from Equations ((TJ) and (J2J). The product of densities J2 Pi(Pi)P2(p2,P3) is determined by a 
similar mixing procedure, as defined earlier, where two like-sign tracks from the same event are 
combined with one track having the same charge from another event with the same multiplicity. 
Finally, the variable Q 3 is calculated from these three tracks. This procedure is similar to that 
given in Reference |2H The ratio (J2 P1P2VP0 is also corrected for detector effects as p 3 /p m i x . 

In our analysis, we use JETSET without BEC and HERWIG to determine C m j x and JETSET 
with and without BEC as well as HERWIG to determine Cdet- These six MC combinations 
serve to estimate systematic uncertainties. The corrections are largest at small Q 3 . At Q 3 = 
0.16 GeV, these corrections to R 3 are C mix ~ 5-30% and Cdet ~ 20-30%, depending on which 
MC is used. These corrections for R 3 and R± 2 are correlated and largely cancel in calculating 

^genuine by Equation (JXJJ) . 

To correct the data for two-pion Coulomb repulsion in calculating P2, each pair of pions is 
weighted by the inverse Gamow factor [22] 

lf exp(270fo) - 1 m n a 
2 ™ = 2^ ' re Vij = ~Q~ ^ ' 

and a is the fine-structure constant. It has been shown |2S] that this Gamow factor is an 
approximation suitable for our purposes. For p 3 , the weight of each triplet is taken as the 
product of the weights of the three pairs within it. For P2P1 we use the same weight but with 
G2{Qij) = 1 when particles i and j come from different events. At the lowest Q3 values under 
consideration, the Coulomb correction is approximately 10%, 3% and 2%, for p 3 , X) pip 2 and 
P2, respectively. 



Results 

The measurements of i? 3 , R\^ and R2 are shown in Figure 121 The full circles correspond to the 
averages of the data points obtained from the six possible MC combinations used to determine 
C mix and Cdet- The error bars, <Ti, include both the statistical uncertainty and the systematic 
uncertainty of the MC modeling, which is taken as the r.m.s. of the values obtained using the 
different MC combinations. This dominant systematic uncertainty is, for Q 3 < 0.8 GeV, about 
a factor 5 to 7 larger than the statistical uncertainty and is correlated between the _R 3 , R\^ and 
R2 distributions of Figure 121 and between bins. Figure shows the existence of three-particle 
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correlations and from Figure [2Jd it is clear that about half is due to two-particle correlations. 
Figure Et shows the two-particle correlations. 

As a check, R 3 , and R2 are also computed for MC models without BEC, both HERWIG 
and JETSET, after detector simulation, reconstruction and selection. For the mixing and 
detector corrections all possible MC combinations, giving non-trivial results, are studied. The 
results of this check are shown in Figure |2] as open circles and, as expected, flat distributions 
around unity are observed. 

Figure shows the genuine three-particle BE correlation function i?f enume . The data points 
show the existence of genuine three-particle BE correlations. The MC systematic uncertainty 
is highly correlated from bin to bin. At Q 3 < 0.8 GeV, it is about a factor 1.5 to 3.5 larger 
than the statistical uncertainty, the higher value corresponding to the lowest Q 3 value used. 
The open circles correspond to MC without BEC and form a flat distribution around unity, as 
expected. 



Gaussian Parametrizations 

A fit from Q 3 = 0.16 to 1.40 GeV using the covariance matrix including both the statistical 
uncertainty and the systematic uncertainty due to the MC modeling, a%, is performed on the 
data points with the commonly used [HIEHEIE] parametrization 

£f nuin °(Q 3 ) = 7 [l + 2A L5 exp(- J R 2 Q2/ 2 )] (1 + sQ 3 ) , (13) 

where 7 is an overall normalization factor, A measures the strength of the correlation, R is a 
measure for the effective source size in space-time and the term (1 + eQ 3 ) takes into account 
possible long-range momentum correlations. The form of this parametrization is a consequence 
of the assumptions that u> — 1 and that \F(Qij)\ = \f\ exp(— R 2 Q^/2), as would be expected 
for a Gaussian source density. The fit results are given in the first column of Table ^ and shown 
as the full line in Figure 

In addition to the MC modeling, we investigate four other sources of systematic uncertainties 
on the fit parameters. Firstly, the influence of a different mixing sample is studied by removing 
the conditions that tracks are replaced by tracks with the same charge and coming from events 
with approximately the same multiplicity. For each of the six MC combinations, the difference 
in the fit results between the two mixing methods is taken as an estimate of the systematic 
uncertainty. The square root of the mean of the squares of these differences is assigned as 
the systematic uncertainty from this source. In the same way, systematic uncertainties related 
to track and event selection and to the choice of the fit range are evaluated. The analysis is 
repeated with stronger and weaker selection criteria, changing the number of events by about 
±11% and the number of tracks by about ±12%. The fit range is varied by removing the first 
point of the fit and varying the end point by ±200 MeV. Finally, we study the influence of 
removing like-sign track pairs with small polar and azimuthal opening angles. The maximum 
deviation that is found by varying the cuts on these angles up to 6°, is taken as the systematic 
uncertainty from this source. The total systematic uncertainty due to these four sources is 
obtained by adding the four uncertainties in quadrature. We refer to this systematic uncertainty 
as a 2 . For all fit parameters, the largest part of the total uncertainty is due to the six possible 
combinations of mixing and detector MC corrections and amounts to 50 to 90%. Table El shows 
the uncertainties for each of the sources for the fit parameters of Equation . 

As a cross-check, the analysis is repeated without the momentum cut of 1 GeV and without 
the cut of 3° on the opening angle of like-sign track pairs. The results agree with those given 



6 



in Table dwell within quoted uncertainties, but the systematic uncertainties are approximately 
twice as large. 

To measure the ratio u, we also need to determine the two-particle BE correlation function 
R 2 {Q)- This is done in the same way as the three-particle BE correlation function. The 
correlation function R 2 is parametrized as a Gaussian: 

R 2 {Q) = 7 [l + Aexp(-i? 2 Q 2 )] (1 + eQ) . (14) 

The parametrization starts at Q = 0.08 GeV, consistent with the study of R3 from Q3 = 
0.16 GeV. The fit result^] are given in the first column of Table El and in Figure Eh- 

If the space-time structure of the pion source is Gaussian and the pion production mechanism 
is completely incoherent, A and R as derived from the fit to Equation (|13jl measure the same 
correlation strength and effective source size as A and R of Equation (fT4"j) . The values of A and 
R are consistent with A and R, as expected for fully incoherent production of pions {uj = 1). 
Using the values of A and R instead of A and R in Equation (|13|). which is justified if u — 1, 
results in the dashed line in Figure EK- It is only slightly different from the result of the fit to 
Equation (fl"3|) . indicating that uj is indeed near unity. 

Another way to see how well ^| enume corresponds to a completely incoherent pion production 
interpretation and a Gaussian source density in space-time, is to compute uj with Equation (jEJ, 
for each bin in Q 3 (from 0.16 to 0.80 GeV), using the measured ^| enume anc [ R 2 derived from 
the parametrization of Equation (jl4j) . The result is shown in Figure 13^- At low Q3, uj appears 
to be higher than unity. 



Extended Gaussian Parametrizations 



However, the assumption of a Gaussian source density is only an approximation, as observed 
in Reference El and confirmed by the x 2 of the fit to Equation ()14j) . Deviations from a Gaussian 
can be studied by expanding in terms of derivatives of the Gaussian, which are related to 
Hermite polynomials. Taking only the lowest-order non-Gaussian term into account, this so- 
called Edgeworth expansion replaces the parametrization of Equation ()14j) by 

(15) 

where k measures the deviation from the Gaussian and H^x) = x 3 — 3x is the third-order 
Hermite polynomial. The fit results for the two-particle BE correlation function with this 
parametrization are given in the second column of Table El 

Using the first-order Edgeworth expansion of the Gaussian, Equation ()15|) . and using Equa- 
tion (jSJ), assuming uj — 1, the parametrization of Equation ()13|) becomes 



R 2 (Q) = 7 [l + Aexp(- J R 2 g 2 )(l + kH 3 {V2RQ)/6)] (1 + eQ) 



Rf~(Q 3 ) =7(1 + 2A L5 exp(- J R 2 Q 2 /2) 
7 I 1 + 2A L5 exp(- J R 2 Q 2 /2) 




1 ~~ r~ ~ Av 
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[l + eQ 3 



H 3 (y/2EQ 3 /2) J 

1 H K 

6 



1.5N 



(16) 



In the second line the approximation is made that Qij = Q3/2. The effect of this approximation 
on ^| enume i s sma ll compared to the statistical uncertainty. The results of a fit to Equation |T6|) 
are given in the second column of Table d The uncertainties are summarized in Table El 

2 - ) Due to the use of a different fit range, these fit results differ from those found in Reference 1241 The same 
fit range gives similar results. 
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For both ^| enmne anc } ^2, a better % 2 /NDF is found using the Edgeworth expansion, and 
the values of A and A are significantly higher, as shown in Tables Q and El and in Figures Eb 
and Et- The values for A and R are still consistent with the corresponding A and R, as would 
be expected for a fully incoherent production mechanism of pions. 

In FigureEb; as in Figure^; we observe good agreement between the fit of f£| enumc using the 
parametrization of Equation (|16j) and the prediction of a completely incoherent pion production 
mechanism, derived from parametrizing R2 with Equation (|15jl . over the full range of Q3. In 
Figure EJd, no deviation from unity is observed for the ratio u. This indicates that the data 
agree with the assumption of fully incoherent pion production. 

Fits to samples generated with JETSET with BE effects modelled by BE or BE 32 E1 JH] 
result in values of R in agreement with the data but in significantly higher values of A. This 
confirms the observation in Figure d that the standard BE implementations of JETSET over- 
estimate the genuine three-particle BEC. 



Summary 

Three-particle, as well as two-particle Bose- Einstein correlations of like-sign charged pions have 
been measured in hadronic Z decay Genuine three-particle BE correlations are observed. The 
correlation functions are better parametrized by an Edgeworth expansion of a Gaussian than 
by a simple Gaussian. Combining the two- and three-particle correlations shows that the data 
are consistent with a fully incoherent production mechanism of pions. 
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parameter 


Eq. © 


Eq. PI) 


7 
A 

rt, im 

e, GeV- 1 
R 

X 2 /NDF 


0.96 ±0.03 ±0.02 
0.47 ±0.07 ±0.03 
U.oo ± U.UO ± U.Uo 
0.02 ±0.02 ±0.02 

29.9/27 


0.95 ±0.03 ±0.02 
0.75 ±0.10 ±0.03 

U. ( z it U.Uo ± U.Uo 

0.02 ±0.02 ±0.02 
0.79 ±0.26 ±0.15 
17.7/26 



Table 1: Values of the fit parameters for the genuine three-particle BE correlation function 
i?| enuine , using the parametrizations of Equations (fT3*|) and (fT6*|) . The first uncertainty corre- 
sponds to Ox, the second to defined in the text. 



parametrization 
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Eq. JIB 
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fit parameter 
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R, fm 


e, GeV" 1 
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A 


R, fm 


e, GeV" 1 


K 


<7i (stat.±modeling) 


0.029 


0.071 


0.056 


0.022 


0.031 


0.103 


0.078 


0.024 


0.26 


mixing 


0.004 


0.006 


0.009 


0.007 


0.010 


0.009 


0.011 


0.011 


0.04 


fit range 


0.008 


0.019 


0.020 


0.013 


0.010 


0.022 


0.017 


0.019 


0.14 


track/event sel. 


0.010 


0.013 


0.012 


0.008 


0.011 


0.016 


0.012 


0.007 


0.10 


5(j) + 56 cut 


0.013 


0.014 


0.012 


0.009 


0.014 


0.017 


0.010 


0.008 


0.11 


0"2 


0.019 


0.028 


0.028 


0.020 


0.023 


0.033 


0.026 


0.024 


0.15 



Table 2: Contribution to the uncertainty on the fit parameters of the parametrizations of 
Equations (fTo^ and ()16j) . respectively. The first uncertainty corresponds to ai, the others 
added in quadrature give cr 2 . 
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Eq. dUD 
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R, fm 

e, GeV' 1 

K 

X 2 /NDF 


0.98 ±0.03 ±0.02 
0.45 ±0.06 ±0.03 
0.65 ±0.03 ±0.03 
0.01 ±0.01 ±0.02 

60.2/29 


0.96 ±0.03 ±0.02 
0.72 ±0.08 ±0.03 
0.74 ±0.06 ±0.02 
0.01 ±0.02 ±0.02 
0.74 ±0.21 ±0.15 
26.0/28 



Table 3: Values of the fit parameters for the two-particle BE correlation function, R 2 , using 
the parametrizations of Equations (fT4^) and (fT5|). The first uncertainty corresponds to o"i, the 
second to 
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Figure 1: Normalized distributions of the sum of the difference in (a) azimuthal and (b) polar 
angle of pairs of tracks in a triplet, I] 50 and ^56, and of (c) Q 3 . Data, JETSET, with 
and without BE effects, and HERWIG are displayed. The ratios between the data and MC 
distributions are shown in (d), (e) and (f). 
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Figure 2: (a) The three-particle BE correlation function, R3, from Equation (b) the 

contribution of two-particle correlations, i?^ = (2 P2P1)/ Po — 2, and (c) i?2 from Equation (JHJ). 
The full circles correspond to the data and the error bars to o\ (see text). The open circles 
correspond to the results from MC models without BEC. In (c) the dashed and full lines show 
the fits of Equations (|T4|) and (JEJ), respectively. 
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Figure 3: The genuine three-particle BE correlation function Rf enmne , Equation (jllj) . The full 
circles correspond to the data and the error bars to o\. The open circles correspond to results 
from MC models without BEC. In (a) the full line shows the fit of Equation (|13|). the dashed 
line the prediction of completely incoherent pion production and a Gaussian source density in 
space-time, derived from parametrizing R 2 with Equation (fTlj) . In (b) Equations (fTBj) and (JT3J) 
are used, respectively. 



15 



2.25 
2 

1.75 
1.5 
1.25 
1 

0.75 
0.5 
0.25 





Gaussian 

j i i i i i i i i 

0.2 0.4 0.6 

Q q [GeV] 



0.8 



2.25 



1.75 
1.5 
1.25 
1 

0.75 
0.5 
0.25 




(b) 



L3 



Edgeworth 

~ i i i i i i i i i 

0.2 0.4 0.6 

Q q [GeV] 



0.8 



Figure 4: The ratio u as a function of Q3 assuming R2 is described (a) by the Gaussian, 
Equation (JEJ), and (b) by the first-order Edgeworth expansion of the Gaussian, Equation (fT5j) . 
The error bars correspond to 0\. For completely incoherent production, u = 1. 



16 



